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NOMENCLATURE 


collecting  disc  area 

ion  flux  coefficient 

di ameter 
electric  field 
charge  on  electron 
mobility 

Boltzmann's  constant 
mass 

reduced  mass 
number  density 
pressure 

collision  cross  section 

speed  ratio 

temperature 

velocity 

thermal  velocity 
drift  velocity 

ro/re 

V*e 

mean  free  path 
Debye  length 
polari zability 


NOMENCLATURE  (CONT.) 


a  particle  diameter 

<f>  electric  potential 


Subscripts 


"  free  steam  conditions 

e  electron 

i  ion 
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1.  INTRODUCTION  AND  SUMMARY 


The  objective  of  contract  number  FI 9628-77-C-0020  has  been  to 
establish  a  connection  between  the  properties  of  the  ambient  ionosphere 
and  the  strength  and  chemical  state  of  the  flux  entering  a  rocket-borne 
mass  spectrometer.  More  specifically,  the  primary  concern  of  this  study 
was  the  description  of  the  collection  of  negative  ions  by  a  flat-faced 
cylinder,  approximately  15  cm.  in  diameter,  travelling  through  a  weakly 
ionized  gas.  The  approach  taken  was  to  use  the  Direct  Simulation 
Monte  Carlo  Method  originated  by  G.  A.  Bird.  This  method  was  developed 
initially  for  neutral  gases  and  has  been  extended  to  a  plasma  in  the 
presence  of  a  charged  body  with  a  consistent  electric  field.  Since  the 
electric  potential  distribution  requires  a  solution  of  Poisson's  equation 
the  charged  particle  calculations  require  an  iterative  procedure  between 
the  Monte  Carlo  results  and  a  finite  difference  solution  of  Poisson's 
equation.  This  approach  was  successfully  developed  and  applied  under 
previous  contracts  to  predict  the  collection  of  positive  ions.  The 
results  of  these  calculations  along  with  the  details  of  the  direct 
simulation  Monte  Carlo  procedure  are  presented  in  Vogenitz  (1973), 

Sugimura  and  Vogenitz  (1973,1975),  and  Sugimura  (1976). 

In  the  present  study,  the  computer  code  developed  for  positive  ion 
collection  has  been  modified  to  predict  the  flux  of  negative  ions  to  a 
rocket-borne  mass  spectrometer.  The  technical  approach  and  formulation 
of  the  problem  is  outlined  in  Section  2.  The  interaction  between  ions 
and  neutral  atoms  is  addressed  in  Section  3.  Section  4  presents  a 
parametric  variation  of  applied  potential,  Debye  length,  and  Knudsen 
number  on  the  flux  of  negative  ions.  Section  5  investigates  the  simulation 
of  ionospheric  chemistry  as  a  collision  process  and  includes  the  effects 
of  internal  degrees  of  freedom.  Finally,  the  conclusions  are  presented 
in  Section  6  where  it  is  noted  that  the  influence  of  angle  of  attack 
was  not  investigated  due  to  the  excessive  computational  requirements 
associated  with  the  3-dimensional  flowfield  and  electric  field. 


2.  TECHNICAL  APPROACH 


2.1  PROBLEM  FORMULATION 

In  the  present  study,  the  negative  ion  flux  to  a  rocket  borne  mass 
spectrometer  is  investigated.  The  collection  of  negative  ions  requires  a 
positive  collecting  potential  on  the  front  face.  This  positive  potential, 
however,  also  attracts  the  electrons  which  will  reduce  the  collecting 
potential  and  the  effectiveness  of  the  instrument.  The  sampling  of 
negative  ions  requires  a  circular  double  disk  configuration;  Sherman  and 
Parker  (1970),  on  the  front  face  to  repel  the  electrons  while  still 
collecting  negative  ions.  The  geometry  and  the  potential  distribution 
along  the  front  face  of  the  instrument  payload  will  be  approximated  as 
shown  in  Figure  1.  Although  the  face  potential,  would  be  greater 
than  the  mask  potential,  $  ,  the  mask  has  greater  area.  Therefore, 
only  the  proper  ratio  of  $  ,  to  will  allow  the  ions  to  pass  through 
the  collecting  orifice  for  a  given  set  of  free  stream  and  boundary  conditions. 

The  most  important  assumption  of  the  Monte  Carlo  simulation  regards 
the  treatment  of  the  electrons.  For  the  positive  ion  collection,  the  front 
face  is  negatively  charged  which  repels  the  electrons  and  attracts  the 
positive  ions.  The  negative  front  face  potential  was  assumed  large 
enough  so  that  the  electrons  could  be  taken  to  be  in  equilibrium  with  the 
local  potential  and  determined  by  the  Boltzmann  distribution 


.  e(e*ATe) 


where 


Ng  =  electron  number  density 
N„  =  free  steam  value 

QOO 

<j>  =  electric  potential 

Tg  =  electron  temperature 

e  =  charge  on  electron 


This  assumption  was  shown  by  Vogenitz  (1973)  to  be  accurate  and  was  then  used 
throughout  the  previous  studies.  The  great  simplification  of  this  assumption 
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was  that  the  electron  trajectories  could  be  ignored  and  only  the  ion 
motion  was  computed  with  the  only  coupling  between  the  ion  and  electron 
distributions  occurring  through  the  electric  potential  determined  from 
Poisson's  equation 


«c  =  positive  collection  potential  for  negative  ion 

<t  =  negative  electron  mask  potential 
e 

Figure  1.  Front  Face  Geometry 

(Ri+-Re>  (2) 

where 

$  =  =  non-dimensional  electric  potential 

*  oo 

Aw  =  mean  free  path 
Arv  =  Debye  length 
+  N'  + 

fJ.  =  =  positive  ion  number  density  ratio 

i« 

Ne 

=  ip  =  electron  number  density  ratio 
e»  3 

MniWuiMfi -r<ir  if^irf  -i*--  •  •  -  • 


For  the  collection  of  negative  ions,  this  equilibrium  assumption  for  the 
electrons  is  not  rigorously  justified.  However,  it  was  necessary  to 
make  this  same  assumption  since  the  process  of  explicitly  simulating 
electron  trajectories  was  found  to  involve  prohibitive  computational 
complexity.  This  complexity  arises  for  three  reasons:  (1)  significant 
additional  memory  is  required  to  store  the  coordinates  and  velocity 
components  of  the  electrons;  (2)  the  high  thermal  speed  of  the  electrons 
necessitates  a  much  smaller  time  step  in  the  simulation  which  results  in 
a  large  increase  in  computation  time;  (3)  the  explicit  simulation  of 
electron  trajectories  allows  for  greater  statistical  fluctuations  of 
local  charge  densities  thereby  decreasing  the  stability  of  the  iterative 
process  for  determining  an  ion  flowfield  with  a  consistent  electric 

field.  It  should  be  noted  that  this  assumption  has  no  effect  on  the 
calculations  for  large  Debye  numbers  since,  in  that  case,  the  electric 

field  is  not  significantly  affected  by  the  ion-electron  distribution. 

For  smaller  Debye  numbers,  the  results  of  the  simulation  appear  to  be 

realistic  and  hence  it  was  concluded  that  the  approximation  did  not 

produce  a  significant  adverse  effect.  _ _ _ _  - 

For  the  collection  of  negative  ions,  two  hypotheses  regarding  ion 

concentrations  are  appropriate.  The  first  of  these  is  the  same  as 
made  in  the  problem  of  positive-ion  collection,  namely,  the  gas  is 
weakly  ionized  and  hence  the  ions  and  electrons  are  trace  species 
compared  to  the  neutral  atoms.  Thus,  the  flow  of  the  charged  species  is 
affected  by  the  neutrals,  but  not  vice  versa.  Additionally,  the  interaction 
of  charged  particles  occurs  only  through  the  electric  field,  not  through 
binary  encounters.  The  second  hypothesis  is  that  the  concentration  of 
negative  ions  is  small  compared  to  that  of  positive  ions.  The  literature 
on  ionospheric  composition  indicates  that  this  is  indeed  the  case  under 
the  great  majority  of  conditions.  The  benefit  of  this  hypothesis  can  be 
seen  by  noting  that  the  Poisson  equation  appropriate  for  the  general 


where  all  the  variables  are  identical  to  those  defined  in  Equation  (2)  with 
addition  of 
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f).j  =  hj1  =  negative  ion  number  density  ratio 
i« 

If  the  concentration  of  negative  ions  is  negligible,  then  the  term 
involving  R..-  can  be  dropped  and  hence  the  calculation  of  the  negative- 
ion  flux  to  the  cylinder  face  reduces  to  a  process  of  allowing  the  negative 
ions  to  percolate  through  the  neutral -molecule  flowfield  with  a  fixed 
electric  field  determined  solely  by  the  positive  ions  and  electrons. 

This  computational  procedure  for  the  negative-ion  sampling  is  outlined 
more  fully  in  the  following  section. 

2.2  COMPUTATIONAL  PROCEDURE 

The  calculation  procedure  for  the  collection  of  negative  ions  is 
shown  in  the  flowchart  of  Figure  2.  It  consists  of  the  following 
steps : 

(1)  Initialize  the  problem  by  specifying  the  geometry  and  flow  condi¬ 
tions,  e.g., 

•  Vehicle  velocity 

•  Mean  free  path  for  neutral  molecules 

•  Debye  lenth 

•  Electric  potential  distribution  on  the  vehicle. 

(2)  Since  the  ions  and  electrons  are  trace  species,  the  neutral 
molecule  flowfield  can  be  determined  without  considering  the  ions 
and  electrons.  The  neutral  flowfield  computed  by  the  Monte  Carlo 
direct  simulation  method  is  used  as  a  scattering  field  for  the 
ion  and  electron  distributions.  This  calculation  is  performed 
only  once  for  a  given  mean  free  path  and  can  be  used  for  all 
Debye  lengths.  In  the  collisionless  limit,  this  calculation 

is  bypassed. 

(3)  At  this  point,  the  ion  and  electron  number  densities  have  not 
been  determined.  Therefore,  the  first  iteration  for  the  electric 
potential  is  found  by  solving  the  Laplace  equation 

V2$  =  0  (4) 


5 


FIGURE  2.  Flow  Chart  for  Negative  Ion  Collection 
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If  the  Debye  length  is  infinite  (i.e.,  *D  =  »),  this  potential 
is  exact  and  no  further  solutions  of  the  Poisson  equation  are 
needed. 

(4)  The  neutrals  are  kept  stationary  while  the  motion  of  the 
positive  ions  (appropriate  to  a  small  time  interval)  is 
computed  using  the  present  electric  potential  distribution. 

(5)  The  new  charge  density  and  corresponding  electric  potential 
is  calculated  for  the  current  distribution  of  electrons 
and  positive  ions. 

(6)  Convergence  of  the  electron  and  positive  ion  distributions  for  a 
consistent  Poisson  solution  is  assessed.  This  is  accomplished 
by  checking  the  charged  particle  distributions  for  the  current 
results  with  the  previous  iteration  to  see  if  steady  state 

has  been  reached.  If  not,  the  iteration  continues  by  returning 
to  step  (4). 

(7)  With  the  neutrals  and  positive  ions  held  stationary,  the  negative 
ions  are  moved  repeatedly  in  the  fixed  electric  field  to  deter¬ 
mine  the  flux  to  the  instrument  face. 

It  should  be  noted  that  step  (4)  requires  an  understanding  of  the  interaction 
between  positive  ions  and  neutrals.  Additionally,  an  understanding  of 
electron-neutral  interactions  is  needed  if  the  electron  motion  were  to  be 
directly  simulated.  This  subject  of  charged  particle-neutral  interaction 
is  discussed  in  the  next  chapter. 


3.  CHARGED  PARTICLE-NEUTRAL  INTERACTION 


Monte  Carlo  calculations  performed  for  the  electron  motion  in  a  neutral  gas 
at  densities  where  collisions  are  important  indicated  that  the  drift  velocities 
were  much  higher  than  accepted  measured  values.  After  confirming  the  correctness  of 
the  numerical  procedure,  the  basic  assumptions  for  the  simulation  were  reexamined. 

In  the  past  simple  hard  spheres  were  used  for  the  charged  particle-neutral 
scattering  model  and  the  collision  cross  section  was  computed  from 


where  ap  *  charged  particle  diameter 
aR  -  neutral  particle  diameter 
For  the  ions  (p  *  1) 


and  for  the  electrons  (p  *  e) 


The  simple  hard  sphere  model  used  for  convenience  in  the  earlier  studies  was  justi¬ 
fied  on  the  basis  that  the  gross  features  of  the  ion  collection  process  were  not 
clearly  understood  and  in  the  collection  of  positive  ions  the  electron  trajectories 
were  not  computed.  However,  the  present  task  of  simulating  the  collection  of  neg¬ 
ative  ions  requires  accurate  calculations  of  the  electron  motion. 


The  proper  description  of  the  charged  particle-neutral  encounter  can  be  rep¬ 
resented  by  an  interaction  potential  of  the  form 


V (r)  ■  Ae'fr  +  ^  Cnr‘n 
n  ■  4 
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where  the  first  term  represents  strong,  short  range  repulsive  forces  due  to 
ordinary  neutral  particle  gas  kinetics  and  the  succeeding  terms  represent  attrac¬ 
tive  forces  which  arise  between  the  charged  particle  and  the  polarized  structure 
of  the  neutral. 

In  the  calculations  presented  here  a  4^  power  attraction  suggested  origi¬ 
nally  by  Langevin  was  used.  The  attractive  force,  due  to  polarization  of  the 
neutral  molecule  is  given  by 


where  y  =  polarizability 

e  =  electronic  charge 

The  current  Monte  Carlo  calculations  treat  the  charged  particle-neutral 
encounter  as  a  binary  collison  involving  both  attractive  and  repulsive 
forces.  The  particle  trajectories  due  to  the  repulsive  force  can  be 
separated  into  two  effects:  (1)  isotropic  scattering  at  an  impact 
parameter  inversely  proportional  to  the  relative  velocity  of  the  collision, 
(2)  isotropic  scattering  by  a  hard  sphere  core  of  fixed  radius.  The 
dominance  of  (1),  (2)  or  the  attractive  force  depends  on  the  value  of 
the  impact  parameter.  For  low  energy  encounters  the  attractive  force 
would  dominate,  whereas  for  high  energy  encounters  strong  repulsive 
forces  would  be  most  important  (see  Figure  3). 


The  collision  diameter  a 

o 


* 

where  m 


m  m 


=  JLJ1 


m  +  m 
P 


n 


for  the  attractive  force  is  found  to  be 


u„  =  relative  velocity  of  the  encounter 
r 

This  diameter  can  be  compared  to  the  standard  hard  sphere  diameter  an 
to  illustrate  the  polarization  effect.  For  example,  oxygen  has  a  polar¬ 
izability 


u  -  1.6  x 


in-24  3 

10  cm 
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HARD  SPHERE  CORE 
r  =  PRESCRIBED  CONSTANT 


and  a  hard  sphere  diameter 


a  «  3  x  10"®  cm 
n 


The  reduced  mass  for  the  ion-neutral  interaction  becomes 


m_m 


m_ 


m 


.  n  1  ^  n  >_  y  _  \ 

i  z. — .  -  •  ~  ■5—  -  m  ) 

m1  +  mn  *  '  n 


At  300°K  the  relative  velocity  can  be  approximated  by 

2 


u  2  1.4  x10s  a; 

r  in  c 


Finally  the  ratio  of  the  diameters  is  found  to  be 
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This  result  implies  that  the  collision  cross  section  for  oxygen  ion- 
neutral  encounters  at  300°K  is  approximately  25  times  the  typical  hard 
sphere  value. 

To  investigate  the  results  of  this  collision  model  a  simple  relaxation 
process  has  been  studied.  Consider  a  weakly  ionized  gas  composed  of 
electrons,  positive  ions  and  neutral  particles  which  are  initially  in 
equilibrium  at  a  single  temperature.  The  positive  ions  are  singly 
charged  and  of  mass  equal  to  the  neutral  particles.  A  prescribed  constant 
electric  field  is  instantaneously  applied  and  subsequently  maintained 
and  the  molecular  motion  of  the  gas  is  computed  in  time  as  a  steady 
state  is  approached.  The  physical  situation  to  which  this  corresponds 
could  be  imagined  to  consist  of  a  gas  contained  between  parallel  electrodes 
of  unlimited  extent  and  great  separation.  Attention  is  focused  upon  a 
thin  slab  of  this  gas  within  which  the  spatial  gradients  can  be  ignored. 
Since  all  charged-charged  interactions  are  contained  in  the  electric 
field,  which  is  known  and  constant  for  this  case,  the  electron  and  ion 
motions  are  independent.  Only  the  encounters  of  the  charged  particles 
with  the  neutral  particles  while  under  the  Influence  of  the  electric 


11 


field  determine  the  behavior  of  the  charged  particles. 

Quantities  of  interest  for  this  study  are  the  following: 

1)  the  steady  state  drift  velocity  per  unit  field,  e.g.  the  mobility 

2)  the  relaxation  times 

3)  the  temperature  and  shape  of  the  distribution  function  in  the 
steady  state 

In  the  ion  collection  process  by  a  rocket  borne  mass  spectrometer  biased 
at  a  potential  $  the  parameter  e4>/kT  is  used  to  characterize  the  strength 
of  the  electric  field.  In  the  present  example  of  a  gas  contained 
between  parallel  electroaes  the  parameter  E/P,  field  strength  over 
neutral  gas  pressure,  is  used.  The  acceleration  of  a  particle  with 
charge  e  and  mass  m  in  an  electric  field  E  is  eE/m.  The  energy  acquired 
by  the  particle  from  the  field  between  collisions  is  approximately 
eE*A  where  A  is  the  mean  free  path.  Since  A  is  inversely  proportional 
to  the  density 

1  E 

E‘A  E  •  —  'v  "p'T 

or 


Therefore  the  ratio  of  the  field  energy  to  thermal  energy  is  proportional 
to  E/P.  This  result  assumes  the  field  strength  is  low  and  that  the  masses 
of  the  ion  and  neutral  particle  do  not  differ  greatly  (thus  excluding 
electrons).  For  arbitrary  field  strength  and  ion  mass,  E/P  is  less 
useful.  However,  most  experimental  data  are  presented  in  terms  of  it. 

Results  of  computations  made  with  this  model  is  shown  in  Figures 
4  through  7,  Here  the  product  of  the  mobility  K  =  V^/E  where  =  drift 

velocity,  E  =  electric  field  and  the  pressure,  p,  in  atmospheres  is  shown 
for  the  electrons  in  a  neutral  gas  with  molecular  weight  29.  The  calculations 
compare  adequately  with  the  experiments  which  have  been  influenced 
by  effects  not  modeled  in  the  simulation  such  as  charge  transfer,  quantum 
mechanical  resonance  effects  and  clustering. 
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Figure  4.  Electron  Mobility 


E/p  =  13.1  volte/ cm  mm  Hg 


Figure  5.  Electron  Velocity  Distribution 


4.  NEGATIVE-ION  COLLECTION  RESULTS 


4.1  COLLISIONLESS,  INFINITE  DEBYE  NUMBER  RESULTS 

A  series  of  calculations  have  been  performed  to  predict  the  collection 
of  negative  ions  to  the  probe  geometry  shown  in  Figure  1  for  the  case  in 
which  both  the  Knudsen  number  and  Debye  number  are  infinite.  The  primary 
purpose  of  these  calculations  was  to  investigate  the  effects  of  the  ratio 
of  positive  ($q)  to  negative  (4>g)  front  face  potential.  The  fact  that  the 
Debye  number  is  infinite  is  not  as  restrictive  as  might  be  imagined 
since  it  was  found  that  (just  as  in  the  case  of  positive-ion  collection) 
the  Laplace  solution  was  valid  for  Debye  numbers  as  small  as  1.0. 

The  Laplace  solution  for  the  double-disk  geometry  was  given  by  Sherman 
and  Parker  (1970)  and  can  be  characterized  by  the  following  parameters: 

ratio  of  disk  radii,  a  =  ro/rg 

ratio  of  potentials,  y  =  <i>o/4>e 

Since  the  attracting  potential,  4q,  is  positive  and  the  electron  mask  is 
negative,  it  is  possible  that  the  potential  along  the  stagnation  line 
can  become  negative  and  repel  the  negative  ions.  Sherman  and  Parker 
show  that  the  critical  value  of  y  is  given  by 

_  V  l  -  a2 

~Yc  “  i  -yr^ 

where  this  variation  is  shown  in  Figure  8.  If  |y|  <  |y  |  there  is  a 
repelling  potential  along  some  portion  of  the  stagnation  line,  while  for 
| y |  >  I Yc l  the  potential  along  the  stagnation  line  is  always  positive. 

Calculations  were  made  for  y  slightly  larger  and  slightly  less  than 
the  critical  value  for  several  values  of  a.  The  purpose  of  this  series 
of  calculations  was  to  investigate  the  sensitivity  of  negative  ion  collection 
to  variations  in  the  potential  ratio.  The  other  parameters  were  as  follows: 

Knudsen  number  —  =  1000  (free  molecule  flow) 

Debye  number  ~  =  1000  (Laplace  solution) 

Speed  ratio  S  =  3.1 

Stagnation  potential  4>o  =  200 
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FIGURE  8.  Critical  Potential  Ratio 


An  example  of  these  calculations  is  given  in  Figure  9.  The  result  was 
that  moderate  changes  in  repelling  potential  have  only  a  small  effect 
on  the  ion  flux.  In  terms  of  electron  shielding,  however,  the  smaller 
values  (in  magnitude)  of  potential  ratio  are  more  desireable.  As  a 
balance  between  these  two  considerations,  it  appeared  reasonable  to 
choose  a  =  0.5,  y  =  “5  (yielding  a  small  repelling  potential)  for  the 
subsequent  calculations  investigating  the  effects  of  finite  Knudsen  and 
Debye  numbers.  The  choice  of  a  =  0.5  (with  a  commensurate  value  of  y) 
is  motivated  by  the  results  shown  in  Figure  lOwhere  it  can  be  seen  that 
this  choice  achieves  close  to  the  maximum  stagnation-point  ion  flux 
yet  still  provides  a  reasonably  sized  electron  mask. 

4.2  FLUX  VARIATION  WITH  KNUDSEN  NUMBER 

The  simulations  for  various  Knudsen  numbers  were  performed  with  the 
values  of  the  other  parameters  set  as  previously  specified.  (It  should 
be  noted  that  in  all  the  results  reported  in  this  work,  the  sidewall 
potential  was  assumed  to  be  zero.)  For  finite  Knudsen  numbers,  the  cross- 
section  for  ion-neutral  and  col  1  is  ions .becomes  relevant.  In  the  present 
results,  the  ion-neutral  and  neutral -neutral  cross-sections  were  assumed  to 
be  equal.  [Although  this  assumption  is  not  correct,  it  serves  to  provide 
Knudsen  number  variation  at  minimum  computational  expense.  If  an  absolute 
flux  is  desired,  then  an  accurate  estimate  should  be  obtained  by  basing 
the  Knudsen  number  on  the  ion-neutral  mean  free  path,  rather  than  the  neutral - 
neutral  mean  free  path.  From  this  point  of  view,  we  have  effectively 
increased  the  neutral -neutral  mean  free  path  to  equal  the  ion-neutral  mean 
free  path  (rather  than  vice  versa)  with  the  resultant  effect  of  reducing 
neutral  flowfield  gradients  which  are  limited  anyway  by  cell  dimensions.] 

The  variation  of  the  stagnation-point  ion  flux  coefficient  as  a  function 
of  Knudsen  number  is  shown  in  Figure  11.  It  shows  that  the  effect  of  collisions 
becomes  significant  starting  at  a  Knudsen  number  of  order  10.,  This  conclusion 
is  shown  more  directly  in  Figure  12  by  plotting  the  percentage  of  the  total  ion 
flux  contributed  by  freestream  ions  (i.e.,  those  ions  that  have  not  been 
affected  by  the  body  either  directly  through  collisions  with  the  body  or 
indirectly  throuqh  collisions  with  other  molecules  that  have  been  affected 
by  the  body). 
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FIG.  9 .  ION  FLUX  DISTRIBUTION  AS  FUNCTION  OF  RADIUS  FOR  a 
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This  behavior  of  the  negative-ion  collection  probe  is  not  unlike 
its  behavior  when  positive  ions  are  collected.  A  direct  comparison  of  the 
stagnaton-pointion  flux  coefficients  is  given  in  Figure  13  for  the 
parameters  indicated.  In  both  cases  the  effect  of  collisions  starts  to 
appear  at  a  Knudsen  number  of  about  10.  Note  that  the  stagnation-point 
flux  for  the  negative  ions  is  substantially  higher  at  large  Knudsen 
numbers  than  that  for  positive  ions  because  of  the  strong  focussing 
effect  of  the  double-disk  configuration  (whereas  for  positive-ion 
collection  the  flux  is  much  more  uniform  across  the  face  of  the  probe). 

At  Knudsen  numbers  less  than  unity,  this  focussing  is  greatly  attenuated 
by  collisions.  An  additional  comparison  between  negative  and  positive 
ion  collection  is  shown  in  Figure  14  where  the  fraction  of  freestream 
ions  in  the  stagnation-point  flux  is  plotted  versus  Knudsen  number. 

4.3  FLUX  VARIATION  WITH  DEBYE  NUMBER 

The  simulations  for  various  Debye  numbers  were  performed  with  a 
Knudsen  number  of  1000  and  the  other  parameters  set  as  before.  The 
results  of  these  calculations  are  summarized  in  Figure  15  where  the 
stagnation-point  flux  of  negative  ions  is  plotted  versus  Debye  number. 

The  variation  with  Debye  number  is  qualitatively  similar  to  the  results 
for  positive-ion  collection  except  that  the  roll-off  in  the  flux  with 
decreasing  Debye  number  appears  to  be  more  gradual  for  negative-ion 
collection.  Note  that  results  could  not  be  obtained  for  Debye  numbers 
less  than  0.1  due  to  instabilities  in  the  Poisson  solution  and  in  fact 
the  result  at  Oe  =  0.1  may  be  in  error. 

One  observation  that  applies  to  the  complete  range  of  Knudsen 
and  Debye  numbers  investigated  is  that  there  is  no  flux  of  positive  ions 
to  the  stagnation  point. 
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Negative  Ions 


5.  SIMULATION  OF  IONOSPHERIC  CHEMISTRY 


5.1  INTRODUCTION 

Collision  theory  is  well  established  as  an  elementary  classical  theory 

for  bimolecular  gas  reactions,  and  the  concepts  are  readily  incorporated  in 

the  direct  simulation  method.  It  is  a  phenomenological  approach  in  that  the 

binary  reaction  rate  is  obtained  as  the  product  of  the  kinetic  theory 

collision  rate  for  collisions  with  energy  in  excess  of  the  activation 

energy  and  the  probability  of  reaction  or  steric  factor  P  .  The  steric 

factor  may  be  regarded  as  the  ratio  of  the  reactive  cross-section  to  the 

total  collision  cross-section  and  is  inferred  from  a  comparison  of  the 

theoretical  prediction  with  the  measured  reaction  rate.  The  simplest 

theory  (for  example*)  employs  a  constant  steric  factor  and  considers  only 

the  translational  collision  energy  in  relation  to  the  activation  energy. 

Refinements  of  the  method  employ  a  collision  energy  dependent  steric  factor 

and  allow  for  the  possible  contribution  of  the  internal  energy  of  the 

molecules  in  the  collision.  Collision  theory  may  also  be  extended  to  ter- 

o 

molecular  reactions  (for  example  )  by  assigning  a  lifetime  to  each  binary 
collision  and  regarding  the  triple  or  ternary  collision  as  a  further  binary 
collision  between  the  pair  of  molecules  in  the  initial  collision  and  a 
third  molecule. 

The  reaction  of  most  obvious  interest  in  an  aerodynamic  context  is 
the  dissociation-recombination  of  oxygen  and  nitrogen.  These  are  bi¬ 
molecular  reactions  with  finite  activation  energies  for  dissociation  and 
termolecular  reactions  for  recombination.  The  direct-simulation  method  was 
applied  to  this  reaction  by  Bird.  These  calculations  employed  a  constant 
steric  factor  and  the  results  for  the  dissociation  rate  were  in  agreement 
with  the  theoretical  prediction.  Also,  the  density  and  temperature  depen¬ 
dence  of  the  equilibrium  state  was  in  agreement  with  the  law  of  mass 
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action  for  this  reaction.  The  gas  model  was,  however,  unrepresentative 
of  a  real  diatomic  gas  in  that  there  was  no  internal  energy.  A  further 
difficulty  was  that  computing  time  considerations  ruled  out  calculations 
at  temperatures  small  compared  with  the  characteristic  dissociation  tem¬ 
perature.  Similar  calculations  which  include  the  effects  of  internal  energy 

c 

modes  have  recently  been  made  by  Larsen. 
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Bird**  developed  a  collision  theory  analysis  of  the  dissociation- 
recombination  reaction  to  produce  an  expression  for  the  equilibrium 
state  that  is  identical  in  form  to  the  law  of  mass  action.  If  a  parti¬ 
cular  diatomic  gas  is  to  be  simulated,  a  comprison  of  the  theoretical 
and  experimental  dissociation  rates  gives  the  appropriate  steric  factor, 
while  a  comparison  of  the  theoretical  and  experimental  values  for  the 
equilibrium  degree  of  dissociation  provides  the  binary  collision  lifetime 
parameter.  The  model  corresponding  to  this  extended  theory  was  combined 
with  Borngakke  and  Larsen's  model ^  for  the  internal  energy  modes  and  was 
used  to  compute  the  dissociation  relaxation  zone  behind  a  strong  shock 

O 

wave  in  nitrogen.  The  ability  to  calculate  the  equilibrium  state  was 
particularly  valuable  in  setting  up  the  boundary  conditions  for  the  shock 
wave  simulation.  The  simulation  results  were  in  good  agreement  with 
corresponding  experimental  results. 

A  difficulty  that  was  encountered  in  the  nitrogen  simulation  was 

that  the  negative  temperature  exponent  in  the  rate  equation  was  below  the 

lower  limit  allowed  in  the  theory.  This  had  to  be  allowed  for  by  a 

temperature,  as  well  as  collision  energy,  dependence  of  the  steric  factor. 

The  following  section  presents  a  further  extension  of  the  bimolecular 

theory  to  include  the  contribution  of  the  internal  energy  modes.  As  well 

as  being  more  general,  it  allows  much  larger  negative  temperature  exponents 

in  the  rate  equations  and  would  have  coped  with  the  accepted  value  for 

nitrogen.  Unlike  the  existing  collision  theory  based  on  the  "effective 
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square  term"  approach  ’  ,  it  is  not  restricted  to  temperatures  that  are 
small  in  comparison  with  the  characteristic  temperature  of  the  reaction. 

A  further  extension  is  that  the  model  now  deals  with  inverse  power  law 
molecules  rather  than  being  restricted  to  hard  sphere  molecules.  The 
wider  applicability  of  the  new  model  is  particularly  useful  when  dealing 
with  cases  which  involve  a  very  large  number  of  reactions. 

5.2  DERIVATION  OF  STERIC  FACTORS  FROM  REACTION  RATES 
A  typical  bimolecular  reaction  may  be  written 

A  +  B  C  +  D  ,  (1) 
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where  A,  B,  C,  and  D  represent  separate  molecular  species.  The  rate 
equation  for  species  A  may  be  written 


-  -w  =  ¥T)  Vb  -  ¥T)  nCnD  •  (2) 

where  n  is  the  number  density  with  the  subscripts  representing  the  species, 
t  is  the  time  and  k^(T)  and  kr(T)  are  the  forward  and  reverse  rate  co¬ 
efficients.  The  temperature  T  dependence  of  the  rate  coefficient  is 
generally  of  the  form 


k(T)  =  a  Tb  exp  (-Ea/kT)  , 


(3) 


where  a  and  b  are  constants,  Ea  is  the  activation  energy,  and  k  is  the 
Boltzmann  constant.  The  problem  is  to  find  an  expression  for  the  steric 
factor,  or  collision  probability,  Pr  as  a  function  of  the  collision  energy 
Ec  that  leads  explicitly  to  Eq.  (3)  in  the  context  of  collision  theory. 

The  collision  energy  Ec  is  equal  to  the  sum  of  the  relative  transla¬ 
tional  energy  in  the  collision  and  the  fraction  internal  energy  of  the  two 
molecules  that  may  contribute  to  the  reaction.  The  translational  energy 
is  equal  to  1/2  m^c^.  where  mr  =  m^An^  +  m^)  is  the  reduced  mass  and 
cr  is  the  relative  speed.  The  internal  energy  is  best  dealt  with  by 
specifying  that  the  number  of  degrees  of  freedom  of  the  two  molecules  that 
may  contribute  energy  to  the  reaction  is  ^  and  This  number  may  range 
from  zero  to  the  number  of  degrees  of  the  molecule,  although  there  are 

restrictions  based  on  angular  momentum  considerations  on  the  number  of 

g 

rotational  degrees  of  freedom  that  can  be  active  in  a  reaction.  It  is 
readily  shown  that,  for  an  equilibrium  gas,  the  distribution  function  of 
Ec/kT  is 
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n  -  3 
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(4) 
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where  n  is  the  exponent  of  inverse  power  law  intermolecular  force  and  r 
is  the  Gamma  function.  Collision  theory  assumes  that  the  forward  rate 
equation  is  given  by 

TT  '  -  VAB  /f(lf)  V1^)  •  (5> 

fa 

kT 

where  v^g  is  the  equilibrium  collision  rate  for  a  species  A  molecule  with 
species  B  molecules,  and  is  the  activation  energy  such  that  the  reaction 
probability  is  zero  for  Ec  less  than  Ea>  Now 

VAB  =  nB  ^  • 

where  Oj  is  the  total  collision  cross-section  and  the  bar  denotes  the  mean 
value.  A  comparison  of  Eqs.  (2)  and  (5)  then  shows  that 

oo 

Kn  -  vv/  f(n)  pr<J  (kt)  ■  (6> 

fa 

kT 


For  hard  sphere  (n  =  °°)  molecules,  the  standard  result  is  that 


aTcr 


(-l  +  -2)2  (f) 


1/2 


where  d  is  the  molecular  diameter  while,  for  general  values  of  n» 


aTcr 


=  W, 
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2  (K/mr)n_1  c 


n-5 

n-1 


(7) 

(8) 


where  is  the  cut-off  value  of  the  dimensionless  impact  parameter  and 
K  is  the  constant  in  the  force  law.  The  parameter  Wgm  is  related  to  the 
lower  cut-off  value  that  must  be  chosen  for  the  deflection  angle  in  the 
collision.  The  value  of  Pr  will  depend  on  this  choice  and  the  generalized 
collision  theory  is  unsatisfactory  when  the  overall  aim  is  merely  the 
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calculation  of  the  steric  factor  P  .  However,  this  is  not  a  difficulty 
when  this  factor  is  used  in  the  direct  simulation  method  since  a  cut-off 
value  must  be  specified  anyway. 

If  a  constant  steric  factor  is  employed  in  Eq.  (6),  it  may  be  taken 
outside  the  integration  which  then  yields 
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For  integer  values  of  the  first  argument,  the  incomplete  Gamma  function 

may  be  written  as  a  power  series  in  E./kT.  This  is  effectively  done  in 

®  1 

the  "activation  in  many  degrees  of  freedom"  theory  and  the  "square  term" 
theory.^  All  but  the  leading  term  of  the  power  series  may  be  discarded 
when  kT  «  Efl  but,  while  this  is  satisfactory  for  the  dissociation- 
recombination  reaction,  it  is  out  of  the  question  for  fast  exothermic 
reactions.  This  problem  may  be  avoided  by  the  introduction  of  a  steric 
factor  proportional  to  some  power  of  E  -E  multiplied  by  the  term 

C  di 


This  term  enables  the  incomplete  Gamma  function  to  be  transformed  to  a 
Gamma  function  that  does  not  involve  the  temperature,  and  leads  to  the 
desired  form  of  Eq.  (3)  for  k(T).  Moreover,  this  modification  has  physical 
as  well  as  mathematical  advantages.  It  means  that  Pr  need  no  longer  be 
constant  for  very  large  values  of  Ec  and  the  rise  from  zero  at  Ec  =  Efl 
is  no  longer  discontinuous.  Similar  forms  for  the  steric  factor  have, 
in  fact,  been  proposed  for  these  reasons. 

The  full  expression  for  the  steric  factor  that  yields  Eq.  (3)  when 
substituted  into  Eq.  (6)  with  Eqs.  (4)  and  (8)  is 
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The  symmetry  factor  e  is  equal  to  1  for  unlike  molecules  and  2  when  both 
collision  partners  are  of  the  same  species.  The  use  of  Eq.  (7)  in  place 
of  Eq.  (8)  gives  the  following  simplified  form  of  Eq.  (9)  for  hard  sphere 
molecules: 


/£ i+^2  3\  3 

The  minimum  value  of  b  is  now  -y — ^ —  +  2)  rat^er  than  the  "  2  allowed 

by  earlier  theory  which  does  not  take  into  account  the  contribution  of 
the  internal  energy  modes.  The  extended  theory  would  have  covered  the 
experimental  values  for  nitrogen  that  caused  difficulties  in  the  shock 

O 

wave  calculations.  This  general  behavior  is  consistent  with  the  more 
restricted  "square  term"  theory,  in  that  increasing  negative  value  of  b 
are  associated  with  a  greater  contribution  from  the  internal  degrees  of 
freedom. 

The  use  of  Eqs.  (9)  or  (10)  to  obtain  the  steric  factor  for  use  in 
a  simulation  program  ensures  that  the  reaction  rates  will  agree  with  those 
from  the  continuum  rate  equation  in  flow  situations  where  near  transla¬ 
tional  and  rotational  equilibrium  prevail  and  both  the  discrete  particle 
and  the  continuum  approaches  are  applicable.  The  discrete  particle 
approach,  through  the  direct-simulation  method,  permits  results  to  be 
obtained  for  chemically  reacting  flow  situations  with  marked  non-equilibrium 
in  the  translational  and  internal  modes.  The  critical  question  is  the 
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degree  of  realism  in  the  modeling  of  the  real  flow  in  this  regime.  It 
should  be  kept  in  mind  that  the  above  equations  for  the  collision  proba¬ 
bility  constitute  phenomenological  models  which  are  not  unique  and  which 
have  been  based  largely  on  considerations  of  mathematical  tractability. 

The  ideal  would  be  to  have  tabulations  of  reaction  cross-sections  as  func¬ 
tions  of  the  energy  states  and  impact  parameters  based  on  quantum 
calculations  and  experiment.  The  progress  in  this  direction  has  been 
reviewed  by  Toennies,^  and  very  little  of  this  information  is  currently 
available  for  reactions  of  engineering  interest.  When  a  classical  approach 
based  primarily  on  continuum  rate  data  must  be  used,  the  above  model 
appears  to  offer  the  most  practical  approach  for  chemically  reacting  flows 
in  the  transition  regime.  Hopefully,  the  physical  basis  is  sufficiently 
sound  for  the  predictions  to  be  at  least  qualitatively  valid. 

In  order  to  consider  the  direct  simulation  method  as  an  alternative 
to  conventional  continuum  calculations,  the  following  conditions  should  be 
assessed. 

The  advantages  are: 

•  The  gas  need  not  be  in  an  equilibrium  or  quasi-equilibrium 
state.  Not  only  may  each  species  have  a  separate  temperature, 
but  the  velocity  distributions  functions  may  vary  drastically 
from  the  Maxwellian  distribution  and  the  internal  modes  may 

be  out  of  equilibrium  with  the  translational  mode. 

•  The  calculations  make  direct  use  of  the  reactive  cross- 
sections  as  functions  of  the  molecular  properties  and  collision 
impact  parameters.  The  fact  that  these  do  not  have  to  be 
averaged  to  give  rate  constants  as  a  function  of  temperature 

is  particularly  advantageous  for  electron-neutral  and  electron- 
ion  reactions  for  which  the  cross-section  versus  energy  curves 
often  have  closely  spaced  peaks  and  valleys.  The  vibrational 
states  may  be  regarded  as  reactions.  Radiative  reactions  may 
also  be  included. 

The  disadvantages  are: 

•  The  number  of  representative  simulated  modules  is  generally 
restricted  to  some  thousands  and  the  total  number  of  simulated 
collisions  to  some  millions.  This  means  that  it  is  impractical 
to  simulate  reactions  which  occur  less  frequently  than  once 
every  thousand  elastic  colli  sons  (i.e.,  steric  factors  less 
than  about  10" 3).  This  rules  out  some  very  slow  bimolecular 
reactions  and  generally  rules  out  termolecular  reactions. 
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•  There  are  sampling  difficulties  associated  with  trace  species, 
although  these  may  be  alleviated,  to  some  extent,  through  the 
use  of  weighting  factors  such  that  the  number  of  real  molecules 
represented  by  each  simulated  molecule  varies  from  species  to 
species. 

•  Most  of  the  currently  available  data  is  presented  as  rate 
coefficients  as  a  function  of  temperature,  rather  than  as 
reactive  cross-sections  as  a  function  of  the  collision  energy. 
This  data  must  be  analytically  converted  to  give  the  same 
answer  in  the  quasi-equilibrium  situations  for  which  the  con¬ 
tinuum  approach  is  valid.  (The  magnitude  of  the  steric  factors 
given  by  this  conversion  process  provides  a  useful  check  on 
the  physical  validity  of  the  rate  data.) 

The  above  considerations  lead  to  the  conclusion  that  direct  molecular 
simulation  is  most  suitable  for  rapidly  varying  transient  phenomena  in 
the  upper  regions  of  the  atmosphere.  It  can  be  used  in  situations  where 
the  conventional  continuum  approach  would  be  invalid. 


Illustrative  Example 

The  direct  simulation  method  has  been  applied  to  the  vibrational 

12 

excitation  of  steam.  The  reason  for  choosing  steam  is  that  it  has  an 
anomalously  rapid  vibrational  relaxation  rate  for  V-T  collisions.  The 
scatter  in  the  data  is  such  that  the  temperature  exponent  is  uncertain 
and  a  constant  value  of  3.5  x  10~17m3s-1per  molecule  has  been  used  for 
the  deactivation  rate.  The  characteristic  temperature  of  the  reaction 
is  2293K  so  that  the  activation  rate  coefficient  given  by  the  product  of 
the  deactivation  coefficient  and  exp  (-2295/T).  The  resulting  steric 
factors  from  Eq.  (10)  are  shown  in  Figure  16. 

A  homogeneous  relaxation  case  was  computed  for  initial  conditions 
of  T  =  650K,  n  =  10  23  nf 3  and  n*  =  0.  Here,  n  and  n*  are  the  number  den¬ 
sities  of  the  ground  state  and  first  vibrational  level,  respectively. 

Eq.  (2)  becomes 

&  (?)  ■  *■* » >o6  b  -  m  -  ?]  an 


and  energy  conservation  requires 


STERIC  FACTOR 


Figure  16.  Steric  factors  for  the  activation  and  deactivation  of  the 
first  vibrational  level  in  steam. 
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T 


650  -  918  (n*/n)  . 


(12) 


Figure  17 compares  the  simulation  results  with  a  numerical  solution  of 
eqs.  (11)  and  (12).  The  two  are  in  good  agreement  even  though  this  is  a 
sufficiently  fast  reaction  for  a  substantial  departure  from  local  transla¬ 
tional  and  rotational  equilibrium  to  be  expected.  The  equilibrium  concen¬ 
tration  ratio  is  also  in  good  agreement  with  the  equilibrium  value  even 
though  the  simulation  procedures  based  on  collision  theory  cannot  be 
expected  to  satisfy  detailed  balancing. 

In  order  to  illustrate  the  effects  of  translational  equilibrium  in  a 
transition  regime  flow,  the  steam  model  was  then  applied  to  the  steady 
flow  through  a  screen  or  grid.  This  flow  is  schematically  represented  in 
Figure  18  which  also  defines  the  flowfield  dimensions  and  the  initial  and 
boundary  conditions.  This  is  a  one-dimensional  flow  in  which  the  grid  is 
a  hypothetical  device  which  reflects  a  specified  percentage  of  the  mole¬ 
cules  passing  through  it.  The  initial  temperature  was  chosen  such  that 
the  initial  equilibrium  excitation  is  negligible  and  the  grid  solidity 
is  sufficiently  large  to  produce  a  temperature  with  a  significant  degree 
of  equilibrium  excitation.  The  temperature  profiles  along  this  flow  ahead 
of  the  grid  shown  in  Figurel9are  very  similar  to  those  for  a  normal  shock 
wave.  There  is  a  temperature  jump  at  the  grid  location  and,  since  the 
downstream  Mach  number  in  this  flow  is  barely  supersonic,  there  is  some 
upstream  influence  from  the  open  end.  The  parallel  temperature  Tx  based 
on  the  velocity  components  in  the  flow  direction  overshoots  the  overall 
temperature  T  by  a  significant  margin.  The  normal  temperature  Tn  based  on 
the  other  components  overshoots  T  by  a  small  margin,  while  the  rotational 
temperature  T^  lags  behind  T. 

Figure  20 compares  the  simulation  results  for  the  density  ratio 
n*/(n  +  n*)  with  those  obtained  from  the  solution  of  eq.  (11)  using  the 
values  of  stream  velocity,  density  and  overall  temperature  given  by  the 
simulation.  The  non-equilibrium  efffects  are  significant  and,  as  would  be 
expected  from  the  temperature  profiles  in  Fig.  19,  they  are  greater  when 
the  collision  energy  Ec  is  obtained  from  the  translational  modes  only. 
Other  calculations  showed  that  the  vibrational  overshoot  is  much  larger 
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Figure  17.  A  comparison  of  the  actual  and  simulation  results  for 
vibrational  relaxation  in  steam. 
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Figure  18.  A  schematic  representation  and  data  definition  of  the  test 
case  for  the  flow  of  steam  through  a  grid. 


38 


degree  of  excitation 


700 


,•00 


3 

*- 

2 

t 

s 


•00 


4001 


•  T 

*  T, 

♦  T„ 
O  T, 


'««» 


• 

♦  ,•*'  ooO°O0*,^*IS«» 

♦  •  o  * 


*• 
+ * 


/  ♦*  < 

«“  oO° 


•0079 


-0-09  -0-029 

DISTANCE  <">) 


0  029 


0-09 


Figure  19.  Temperature  profiles  for  the  flow  of  steam  through  a  grid. 
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Figure  20.  Vibrational  excitation  profiles  for  the  flow  of  steam  through 
a  grid. 


39 


if  the  rate  coefficient  is  proportional  to  a  finite  positive  power  of 
temperature,  but  is  not  greatly  affected  by  the  replacement  of  the  basic 
hard  sphere  model  by  an  inverse  ninth  power  law  model. 

The  above  example  illustrates  the  scope  of  the  direct-simulation 
method  when  chemical  data  is  available  only  in  the  form  of  bulk  reaction 
rates.  Approximate  quantum  cross-section  predictions  are  available  for 
T-V  reactions  such  as  that  for  steam  in  the  second  example.  The  direct- 
simulation  approach  would  also  be  particularly  useful  for  the  reverse  task 
of  deducing  rate  coefficients  or  even  flowfields  from  these  predictions 
in  order  to  make  comparisons  with  experiment. 


5.2  NEGATIVE  ION  COMPOSITION  IN  THE  IONOSPHERE 


A  set  of  chemical  reactions  representative  of  the  important  processes 
in  the  ionosphere  has  been  investigated  by  the  direct  simulation  Monte 
Carlo  method  described  earlier.  In  this  calculation,  a  homogeneous  gas 
mixture  initially  containing  electrons,  atomic  oxygen,  molecular  oxygen, 
and  molecular  nitrogen.  Equation  (10)  is  used  to  convert  the  temperature 
dependent  rate  coefficients  given  below  to  collision  energy  dependent 
steric  factors. 

k(cm3/s) 


Formation:  e  +  02  0"  +  0 

e  +  02  +  02  Og  +  02 

Exchange:  0“  +  Og  -*■  0  +  0£ 

Og  +  0  -*■  o2  +  0“ 

Detachment:  0"  +  0  -*•  02  +  e 

©2  +  ©2  0  +  Og  +  e 


<  10"16  for  T  <  2000°K 

(1.4*2)  xlO-»(ip)«p (^)(4-) 
for  (195  <  T  <  600) 


<  10' 


■10 


av  _  in 
<  10  IU 


k(cm3/s) 


1.7+  0.3x10' 


•10 


2.3  x 10”9  exp 


Hf) 


for  T .  <  2000 
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°2  +  0  °3  +  e 


3.2  x 10'10 

V 

®2  +  ^2  +  ®2  +  e  (2.7  +  .3)  x  10  (300)  exP  (  T  ) 

for  (375  <  T  <  600) 

-1?  /  T  / -9990 

^2  +  ^2  **  ®2  +  ^2  +  e  (1.94  +  .4)  x 10  (300)  exP  (  7 

for  (375  <T<600) 

0"  +  N2  0  +  N2  +  e  3.4  x  10"9  exp 

Most  of  the  rate  coefficients  were  obtained  from  the  DNA  Reaction  Rate 
Handbook  (DNA. 19484),  March  1972. 

Two  formation  reactions  are  considered.  One  of  these 

e  +  02  +  02  -*>  O2  +  02 

is  a  termolecular  reaction  with  a  rate  coefficient  that  indicates  that 
it  occurs  only  once  every  107  binary  02>  02  collisions  at  a  number  density 
of  109  to  1010  cm-3.  The  other  reaction 


e  +  02  -*■  O'  +  0 

is  bimolecular  but,  being  dissociative,  has  a  high  activation  energy  of 
approximately  3.6  eV.  Above  this  energy,  the  steric  factor  is  of  the 
order  of  0.2.  The  exchange  and  detachment  reactions  have  steric  factors, 

_3 

in  the  energy  range  of  interest,  that  range  from  the  order  of  10  to 
unity. 

The  test  case  that  has  been  chosen  is  an  atmosphere  with  a  number 
of  10“  cm*3  (approx.  200  km),  with  initial  mole  fractions  of  0.08  electrons, 
0.74  atomic  oxygen,  0.14  molecular  oxygen  and  0.04  molecular  nitrogen, 
and  at  a  temperature  of  2000° K.  The  energy  of  the  electrons  is  instan¬ 
taneously  raised  to  3.6  eV  and  the  reactions  are  followed  for  0.4  seconds. 
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The  composition  is  shown  in  Figure  6  as  a  function  of  time.  The 
mole  fraction  of  0  increases  to  almost  0.04  around  0.12  seconds  and  then 
decreases.  The  fraction  of  0 ^  continues  to  increase  over  the  period  of 
the  calculation.  Figure  7  shows  the  species  temperature  over  the  same 
period.  They  are  still  a  way  from  equilibrium  at  the  end  of  the  calculation 
(but  this  is  a  comparatively  short  calculation).  In  Figure  8  the 
differences  in  the  translational  and  internal  temperatures  shown  as  a 
function  of  time  indicate  the  nonequilibrium  effects  of  the  internal 
degrees  of  freedom. 
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6.  CONCLUSIONS 


The  efforts  to  extend  the  direct  simulation  Monte  Carlo  method  to 
the  collection  of  negative  ions  have  been  shown  to  be  more  difficult 
than  anticipated.  As  a  result,  it  has  been  necessary  to  introduce  a 
simplifying  assumption  for  the  behavior  of  the  electrons  at  finite 
Debye  numbers.  Despite  these  difficulties,  some  useful  results  have  been 
obtained.  In  particular,  it  has  been  found  that  realistic  probe  potentials 
(e.g.,  =  200,  $e  =  -40,  a  =  0.5)  can  produce  a  high  stagnation-point 

flux  of  negative  ions  while  completely  rejecting  positive  ions.  It 
was  also  discovered  that  this  flux  of  negative  ions  is  strongly  attenuated 
by  collisions  at  finite  Knudsen  numbers. 

The  application  of  the  direct  simulation  Monte  Carlo  method  to 
describe  the  chemistry  in  the  ionosphere  has  been  shown  to  be  potentially 
a  very  useful  computational  tool.  Although  the  calculations  presented 
were  restricted  to  a  reaction  cell  simulation  where  there  are  no 
spatial  gradients,  the  method  does  describe  the  effects  of  the  electrons 
and  ions  and  nonequilibrium  influences  which  are  not  accounted  for  in 
the  continuum  descriptions.  The  procedure  can  be  used  to  parametrically 
study  the  effect  of  high  energy  electrons  on  the  chemistry  of  the 
ambient  ionosphere  or  as  a  screening  procedure  for  possible  reaction 
chains.  However,  computational  complexities  prevented  the  incorporation 
of  the  chemistry  simulation  into  the  present  ion  collection  code,  although 
this  remains  as  a  possible  activity  for  future  work. 
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